Ground-state reference systems for expanding correlated fermions in one dimension 



00 
O 

o 



oo 



F. Heidrich-Meisner, 1 ' 2 M. Rigol, 3,4 A. Muramatsu, 5 A. E. Fciguin, 6 and E. Dagotto 2 

1 Institut fur Theoretische Physik C, RWTH Aachen University, 52056 Aachen, Germany 
2 Materials Science and Technology Division, Oak Ridge National Laboratory, Oak Ridge, Tennessee 37831, USA ana 
Department of Physics and Astronomy, University of Tennessee, Knoxville, Tennessee 37996, USA 
1 Department of Physics, University of California, Santa Cruz, California 95064, USA 
^Department of Physics, Georgetown University, Washington, District of Columbia 20057, USA 
J Institut fur Theoretische Physik III, Universitat Stuttgart, 70550 Stuttgart, Germany 
6 Microsoft Project Q, University of California, Santa Barbara, California 93106, USA 

(Dated: May 19, 2008) 

We study the sudden expansion of strongly correlated fermions in a one-dimensional lattice, 
utilizing the time-dependent density-matrix renormalization group method. Our focus is on the 
behavior of experimental observables such as the density, the momentum distribution function, and 
the density and spin structure factors. As our main result, we show that correlations in the transient 
regime can be accurately described by equilibrium reference systems. In addition, we find that the 
expansion from a Mott insulator produces distinctive peaks in the momentum distribution function 
at k ~ ±tt/2, accompanied by the onset of power-law correlations. 
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I. INTRODUCTION 

The noncquilibrium properties of strongly correlated 
electron systems are a challenging subject in need of a 
better understanding. While experimental studies in this 
area face difficulties in the solid state context, ultracold 
quantum gases provide a controlled way to address this 
difficult issue. For this reason, recent experiments em- 
ploying out-of-equilibrium cold atom gases in optical lat- 
tices, which allows the realization of model Hamiltonians 
for strongly correlated particles (for a review, see, e.g., 
Ref. [1]), have attracted considerable attention [2-4]. 

Among the fundamental questions recently addressed 
in these experiments is the issue of thermalization in 
isolated quantum systems [2, 5-14]. In the transient 
regime, quantum quenches have been shown to induce 
a collapse and revival of coherence properties [3], and 
transport measurements in different lattice systems have 
unveiled the intriguing consequences of strong correla- 
tions [4] . The important effects of interactions have been 
observed in the expansion of bosons in one-dimensional 
(ID) lattices as well [15-19]. In the expansion from a 
Mott insulator (MI) state, quasi-condensates at finite 
momenta emerge [15], while in the hard-core regime, the 
expansion from a superfluid state leads to the dynami- 
cal fermionization of the bosonic momentum distribution 
function (MDF) [16]. The latter is a generic feature of 
the expansion of harmonically trapped hard-core bosons 
in the absence of a lattice [20]. In addition, it has been 
shown in Ref. [21] that, independently of the initial inter- 
action strength, a freely expanding Lieb-Liniger gas al- 
ways enters a strongly correlated (hard-core like) regime. 
The expansion dynamics of strongly correlated fermions, 
which due to the spin-degree of freedom is expected to be 
richer, has not yet been addressed, and it is the objective 
of this work. 

Concretely, we study the expansion of two-component 
interacting fermions in a ID lattice. The ground state 
physics of these systems is characterized by a Tomonaga- 



Luttinger (TL) state with power-law decaying correla- 
tions at any incommensurate filling. At half-filling, a 
charge gap opens and the system exhibits quasi long- 
range antiferromagnetic correlations [22]. Here, we wish 
to elucidate how the initial state of the system, being ei- 
ther MI or TL, affects the expansion process. Identifying 
distinctive features for the MI is of much interest to ex- 
perimentalists in the search for the fcrmionic MI state. 
However, our main objective is to shed light on the re- 
lation, if any, between these out-of-equilibrium systems 
and their equilibrium counterparts. As the main result 
of this work, we provide evidence that correlations mea- 
sured in noncquilibrium are quantitatively described by 
appropriately chosen equilibrium reference systems. 

The outline of the paper is the following. First, we 
describe the model and the numerical procedure in Sec. 
II. Section III contains our results on the time-evolution 
of density profiles and the momentum distribution func- 
tion for both TL and MI initial states. In Sec. IV, we 
investigate what the possible relation to equilibrium sys- 
tem is, and we present a comparative analysis of spin 
and charge correlation functions. We also comment on 
the validity of our findings in other models, such as the 
Hubbard chain with a nearest-neighbor repulsion, which 
renders the model nonintcgrable. We conclude with a 
summary of our results contained in Sec. V. 



II. MODEL AND NUMERICAL METHOD 

The noncquilibrium dynamics is analyzed using the 
adaptive time-dependent density-matrix renormalization 
group method (tDMRG) [23]. We consider the ID Hub- 
bard model with nearest-neighbor hopping t and an on- 
site Coulomb repulsion U : 
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FIG. 1: Expansion from a: (a) MI (ni n it=l; main panel: 
MDF; inset: density); (b) TL (n ini t=0.6); both at U = 8t 
and plotted at times r = 0, 4, 8. Note that at r > 8 and 
before the right boundary is reached, the MDF exhibits only 
small changes. Dotted line in (b): MDF of a reference system 
(see text in Sec. IV A for details) with {ni) Te i, T = (^;( r = §))• 
Inset in (b): Decay of one-particle correlations during the 
expansion of the MI. The t~0 curve is for a N=50 system 
at half-filling. Dotted vertical lines in (a) and (b) denote 
k = it/2. 



c l a ( c l,cr) is a fermion creation(annihilation) operator act- 
ing on site I, with (pseudo-)spin index a =1,4, n; )CT = 
Cj cr c i,<r is the corresponding density operator, and we de- 
fine ni = n; jCr . N denotes the number of sites, a is the 
lattice constant, and open boundary conditions are im- 
posed. We prepare an initial state with a filling ni n j t that 
is non- vanishing in only a portion of the system by apply- 
ing a confining box-potential H con f = XwLi e i n i- Hence, 
we have H = Hq + H con {, with e; = 10 6 t for I $ [lo,h] 
and ei = otherwise. At time r = 0, we turn off H con f. 
t is given in units of 1/t, we set H to unity. 

In our tDMRG runs we use a third-order Trotter- 
Suzuki time-evolution scheme with a time step of At = 
0.005. The discarded weight during the time-evolution 
is kept below 10 -8 . To simulate the longest time scales 
possible on a given system size before the particles are re- 
flected at the boundaries, we select an asymmetric set-up 
and, hence, particles can only expand into one direction. 
We have checked that the same overall picture is observed 
in symmetric set-ups (see also Ref. [15]). 



III. THE MOMENTUM DISTRIBUTION 
FUNCTION 

We first discuss the properties of the MDF rife, com- 
puted from n k = (1/N) J2i, m .a- ex P[~^Q~ m )]Phn^ where 
Pirn = ( c l a c m.a) ls the one-particle density matrix. In 
Fig. 1(a), we show the evolution of nk (main panel) and 
the density (m) (inset) for an initial MI state. The main 
panel reveals a peculiar behavior of rik' as the Mott in- 
sulator melts, a peak develops at a finite momentum k p . 
We further find that, for U larger than the band- width 
W = At, k p closely approaches 7r/2. This behavior resem- 
bles that of hard- and soft-core bosons [15]. Qualitatively, 
we understand this in terms of an energy argument: in 
the MI state with U 3> W, the total kinetic energy is 
close to zero. Hence, particles emitted into an empty 
lattice have a small average kinetic energy corresponding 
to a momentum 7r/2. As the Fermi statistics prohibits 
quasi-condensation into a single momentum state, nk be- 
comes a broad function around k p w 7f/2. 

While the initial MI state is characterized by an 
exponential decay of one-particle correlations, i.e., 
l/ 9 taH ex P(~l^ ~ 771 I/O; C=const, we find that during the 
expansion, the system develops power-law correlations. 
In the inset of Fig. 1(b), we compare the \pf m \ of a MI in 
equilibrium with the correlations that emerge during its 
expansion, measured within the moving cloud. The inset 
reveals the weak decay of correlations during the expan- 
sion, consistent with a power law. One may associate 
the dynamical emergence of this power law with a met- 
allization of the moving cloud, which, after the melting 
of the MI, starts behaving as an inhomogeneous metal. 
As of now, our numerical analysis is restricted to a small 
number of particles and time scales of r ~ 15 only, which 
prevents us from extracting, e.g., exponents of the power 
laws. Note, though, that in the case of free fermions 
expanding from an insulating state with ni n j t = 1, i.e., 
a state with no off-diagonal correlations, the emergence 
of power laws has been established for a large number 
of particles and hence over substantially larger distances 
than in the present work [24]. 

In the main panel of Fig. 1(b), we show the evolution 
of the MDF starting from a TL state with n- m i t < 1. 
In this case, the initial state has a well-defined Fermi 
momentum and a power-law decay of correlations [22]. 
Such decay is preserved during the expansion. Moreover, 
rik also exhibits a peak, but at a momentum k p < tt/2 (k p 
increases as n- mi t — > 1). Another property of this peak, 
distinguishing it from the peak formed after the melting 
of the MI, is that it exhibits a much sharper edge at the 
large momentum side, reminiscent of a Fermi edge. 

From the previous analysis, we conclude that if rik 
could be experimentally studied during the expansion in 
the strongly correlated regime, then the emergence of 
peaks at k = ±ir/2 in the fermionic MDF would serve 
to identify the presence of a Mott insulator in the initial 
state. The experimental challenge is to independently 
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FIG. 2: (Color online) Time evolution of (a) the kinetic energy T^ n and (b) the potential energy -Emit for the expansion from 
an initial TL state with U = St, n-mit = 0.6 (dashed lines) and the corresponding reference systems (squares). Panel (a) further 
shows the kinetic energy of the cloud and of the reference systems as estimated from Eq. (5) (thick dotted line and stars, 
respectively) as well as our approximate results for the kinetic energy T^ 1 associated with the center-of-mass motion of the 
cloud [solid line in (a)]. See the text in Sec. IV C for definitions and details. 



control the trapping potential and the lattice [25-27]. 
This has been achieved in the experimental study of dis- 
ordered ultracold Bose gases in both ID optical lattices 
[26] and homogeneous ID systems [27]. 

At this point we would like to emphasize that the 
physics of our expanding system is different from the one 
found in theoretical studies of strongly correlated systems 
in ID lattices undergoing a relaxation following a quan- 
tum quench [5, 9, 10]. In the latter, correlations have 
been found to decay faster than with a power law (some- 
times clearly exponentially) [5, 7, 9, 10], while in our 
moving clouds we find power-law decaying correlations 
[inset of Fig. 1(b)]. In addition, after the relaxation to a 
steady state, one can ask the question of what statisti- 
cal ensemble may best describe physical observables, but 
here we are solely concerned with the transient regime, 
i.e., a regime in which statistical ensembles do not pro- 
vide us with insights into the behavior of physical ob- 
servables. 



IV. GROUND-STATE REFERENCE SYSTEMS 

A. Construction of reference systems 

Our results thus far have singled out a noticeable prop- 
erty of these systems during their expansion: indepen- 
dently of the initial ground state, power laws are ob- 
served in the noncquilibrium dynamics. In ID systems 



in equilibrium, power-law correlations are only seen in 
the ground state, as finite temperatures introduce a cut- 
off at large distances, followed by an exponential decay 
[22]. Hence, one may wonder whether a system out of 
equilibrium can in some way resemble the ground state 
of a system in equilibrium. A natural choice for such a 
reference state is the ground state of a system that has 
exactly the same density distribution as the time evolving 
state, i.e., 

(ni)rei,T = («i( T )) ) ( 2 ) 

for all sites /. Hence, a reference system has to be deter- 
mined at each given time. We construct such reference 
states by self-consistently computing a set of onsite en- 
ergies ei in 

N 

Href = Hq + ^ e l n l ( 3 ) 
1=1 

with Ho from Eq. (1) such that at a desired time, the 
density profile (hj(t)} is reproduced, while keeping t and 
U fixed. Once the density has converged within an error 
of 

5n T = \(ni(r)) - (ni} r ef, T |/X)M r )} < 10_3 - ( 4 ) 
i i 

we compare quantities of interest in both systems. 



4 



B. Spin and charge structure factor 

We now turn to the comparative analysis of corre- 
lation functions. We compute the spin-spin (Sk) and 
density-density (Nk) structure factors, which are the 
Fourier transforms of the spin-spin (Si, n = (SfS^)) 
and density-density (iV; m = (riin m ) — (ni)(n m )) corre- 
lations, respectively. The spin operator is defined as 
Sf = (njj — n;,jJ/2. In equilibrium and for a homo- 
geneous system, Sk peaks at 2kp while Nk exhibits a 
kink at 4fcp [22], where kp = irn/2 is the Fermi mo- 
mentum. Consistently, for the two cases n- m i t = 1 and 
"init = 0.6 shown in Figs. 2(a) and 2(c), respectively, 
Sk(T = 0) [dotted lines] peaks at k = n and k = n0.6, 
while Nk in the case of n- m it = 0.6 has a weak kink at 
k = 27r0.6 [dotted line in Fig. 2(d)]. During the expan- 
sion, the peak in Sk shifts to smaller momenta and the 
maximum is less sharp, as shown in Figs. 2(a) and 2(c) 
[solid lines]. Qualitatively, we understand this behavior 
in terms of the decrease of the average density during the 
expansion into the initially empty lattice, giving rise to a 
shift of the 2kp peak in Sk and a broadening due to the 
inhomogeneity. Further, we propose an operational defi- 
nition of a Fermi-momentum k F in the expanding clouds 
by taking the position of the peak in Sk, yielding 2k T F . 
This supports the use of the term "metallization" for the 
process that fermions escaping from a MI undergo. 

The density correlation Nk docs not show any par- 
ticular features during the time-evolution, as the kink 
is washed out due to the inhomogeneity. Nk increases 
monotonously with time, reflecting an increase of the 
overall charge fluctuations, due to the closing of the 
charge gap as the MI melts. 

Our most remarkable finding, and thus the key result 
of this work, is the excellent agreement seen in Fig. 2 
for Sk and Nk between the expanding cloud - a genuine 
nonequilibrium situation - and the inhomogeneous ref- 
erence systems, which are in their ground state (circles 
in Fig. 2). We are therefore led to conclude that, dur- 
ing the expansion, spin and charge correlations out of 
equilibrium are, to a very good approximation, the same 
functionals of the density as the ones in equilibrium sys- 
tems in the ground state. 



C. Time-evolution of kinetic and potential energy 

Intuitively, we expect that our way of preparing the 
reference systems must yield properties similar to those 
of the moving clouds on short time scales. However, the 
agreement between the clouds and the equilibrium sys- 
tems exists both at short and long times, and is thus 
preserved during the expansion into the empty lattice. 
A noticeable difference exists between nk of the moving 
clouds and nk of the corresponding reference systems. As 
discussed before, nk of the expanding cloud has a finite- 
momentum maximum, which is also present in the sym- 
metric expansion, while nk of the reference system in its 



-6- 



(a) U=8t, n. =0.6 (TL) 



Center-of-mass 
motion of cloud 




time-dep. 

□ ref. system 



^kin ^k^k 



-o □ 

12 



FIG. 3: (Color online) Time evolution of (a) the kinetic energy 
Tki n and (b) the potential energy Emit for the expansion from 
an initial TL state with U = 8t, n-mit — 0.6 (dashed lines) and 
the corresponding reference systems (squares). The plot fur- 
ther shows the kinetic energy of the cloud as estimated from 
Eq. (5) [thick dotted line and stars] as well as our approxi- 
mate results for the kinetic energy T^n' associated with the 
center-of-mass motion of the cloud (solid line). See text in 
Sec. IV C for definitions and details. 



ground state is symmetric around k = [sec dotted line 
in Fig. 1(b)]. Hence, observables related to rik cannot be 
accounted for with this procedure. 

Our previous findings on the spin and density structure 
factors may seem puzzling: at any given time, the sum 
of the kinetic, Xki n , and interaction energy _Ei nt of the 
time-evolving state is much higher than the energy of the 
reference system, which, having the same density profile, 
is in its ground state. This difference grows with time. In 
the case of ninit =0.6 and U = 8t, we have Eq = T^ n + 
E- mt = — 6.86t, which of course is a constant in time. At 



t = 8, this splits into T^m = —t Yli 
-6.97i and E iat = Uj2i( n i,1 n i.i) ; 



h.c. 



0.1 It. In contrast, 

'ref 



for the reference system, we find T^f* = — 10.93£, and 
£^ t f = 0.18t, adding up to E r ci = -10.75*. Thus, the 
main difference is due to the kinetic energies. 

We argue that both systems can be related by a 
Galilean transformation and, thus, understand why their 
structure factors are similar. This means that the dif- 
ference between the kinetic energies T kin and T£f[ is 
mainly due to the average momentum of the moving 
cloud [ko = E fe W£ fc n k ], i.e., that « Tg*, the 
latter being the kinetic energy of the particles in a ref- 
erence frame moving with the cloud. To prove this, we 
first notice that, using the MDF, the kinetic energy can 
be estimated as 
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k 



(5) 



where e° 



-2t cos k is the dispersion relation in the non- 
interacting case. This assumption leads to T^n sa — 6.97i 
and T^fl w — 10.93t, as estimates for the kinetic en- 
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FIG. 4: (Color online) Relative deviation of the density struc- 
ture factor vs. time [see Eq. (6)], for: (a) ni n i t =0.4, 0.6 and 
a harmonic trap (8 fermions, Vtra P =0.02t, U=8t) and (b) 
nmit = l (U/t—8, 20, oo). We estimate our calculation of SNk, T 
to have an accuracy of ±0.01 in the worst case. Inset in (b): 
density in the leftmost site vs. time for U = 20t. 
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ergy of the expanding and reference systems at r = 8 
("■init = 0.6, U — 8t), respectively. Both values are very 
close to the exact results presented before. We then com- 
pute T k <? n al = Ek n k4-k = -10-15* « T£l at r = 8, 
which corroborates our interpretation: the energy differ- 
ence is mostly due to the finite momentum of the cloud, 
and not due to contributions of the internal kinetic or 
interaction energy. This picture is further supported by 
Fig. 3 that contains the time-evolution for T kin , E int , TJ^ , 
and for the parameters discussed in this section. 



D. Breakdown of the ground-state reference 
system description 



FIG. 5: (Color online) U — oo, 500 particles, ntait = 1. (a) 
Decay of \Nij\ at time r = 600 for j = 20 (solid line with 
squares). The dashed lines are fits to the envelope of \Nij\ 
using f(\i — j\) = ot\i — ■ Note that two regions appear 
that are described by the same exponent /3, but a different 
prefactor a. (b) The same as in (a), but for j — 740 (solid 
line). Note that we plot correlations for i > j only, (c) 
Density-density correlations measured in the reference system 
for r = 600, with j = 20 (solid line with squares) and j = 740 
(solid lines). The latter curve has been offset for clarity. The 
dashed-dotted lines are the fits from panels (a) and (b). (d) 
Density profile at time r = 600. Vertical, dashed lines in (a) 
and (d) mark i = 705, separating the regions with different 
prefactors in the power-law decay of \Nij | at the time r = 600 
considered in this plot. 



It is next important to identify conditions for a break- 
down of the reference-system description. To this end, we 
study the relative difference between the time-dependent 
and the reference systems' density structure factors, 

5N k , T = ]T \N k (r) - JVj^|/X; N k (r)- (6) 

k k 

The corresponding errors in S k are smaller than those in 
Nk, and we thus concentrate on the latter. Let us start 
with the initial TL. We consider two cases: first, nj n it < 1 
in a box trap [see Fig. 2(d)]. Second, as such a set-up 
is more realistic to account for experiments, we follow 
the evolution of fermions escaping from a harmonic trap 
Vtrap X);(0 2n i- From the results displayed in Fig. 4(a), 
our key observation is that 6Nk. T ^ 0.02 remains very 
small in both cases. Hence, for an initial TL state and for 
both iVfc and S k , the description given by the equilibrium 
systems is very good up to the largest times simulated. 

We next turn to the case of an initial MI region, and 
present results in Fig. 4(b) for U = St, 20i, oo. The 
U = oo case is treated with exact diagonalization, after 
mapping the charge sector of our two-component fermion 
system to spinlcss fermions [22]. A behavior similar to 
the TL case is found at times r < 5, with SN kyT ^ 0.04. 
However, for times after the melting of the MI region 



[ T ^ 5, see the inset in Fig. 4(b)], a substantial increase 
of SN k , T becomes evident, as shown in Fig. 4(b). Thus, 
for the MI expansion, reference systems work well only 
up to the point at which the Mott insulator totally melts. 

This deviation of the time-dependent data from the 
reference systems is associated to the appearance of par- 
ticular coherence properties in the portion left behind by 
the moving cloud after the melting of the MI, a feature 
that is not captured by the reference systems. To sub- 
stantiate this interpretation, we present results for the 
decay of Nij in the U = oo limit in Fig. 5. In that limit, 
we are able to exactly treat arbitrary time scales for both 
large number of particles and a large system size. The 
data displayed in Fig. 5 are for 500 spinless fermions ex- 
panding from a MI region into a lattice of 2000 sites. 

In Fig. 5, we measure the density-density correlations 
at a time r = 600, at which the initial Fock state has 
already completely melted. Figure 5(d) shows the den- 
sity profile at this time and from that plot, we see that 
at time r = 600, i ps 705 separates the front of fast 
particles from slower ones, as indicated by the vertical 
dashed line in the figure. Panels (a) and (b) show \Nij\ 
for i > j, and j = 20 and j = 740, respectively, i.e., 
measured from behind and from inside the moving front 
[Fig. 5(d)]. While in the latter case, clearly a unique 
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power-law decay of |-/Vy| is observed, the former case is 
more involved. There, the envelope of the correlator also 
decays according to |iVij| = a\i — yet for i < 705 
and i > 705, a different prefactor a is found. These are 
the two regions in Fig. 5(a) indicated by the arrows. The 
exponent (3 = 2K is universal and expected to be (3 = 2 
since the Luttinger parameter K of spinlcss fcrmions is 
K = 1, but the prefactor - in the ground-state - is es- 
sentially a function of the average density, or the Fermi 
momentum, respectively [22]. The density-density cor- 
relations therefore exhibit a distinctly different behavior 
comparing the moving front of fast particles (i > 705) 
and those left behind (i < 705). It is exactly this step- 
like feature, i.e., the sudden change in the prefactor of 
the power law followed by \Nij\, that is not captured by 
the reference systems. This is revealed in Fig. 5(c), which 
shows the [iVy | as measured in the reference system con- 
structed for time r = 600, for both j = 20 and j = 740. 
The plot includes the fits to \Nij\ from panels (a) and 
(b) [dot-dashed lines]. While for j = 740, \N l:j \ is well 
described by the reference systems, \N^° l \ with j = 20 
does not show the step like feature observed in the mov- 
ing cloud, as the reference systems fail to account for the 
separation of particles moving at different velocities. 

For this reasoning to apply, it is important to realize 
that such a separation of velocities as reflected in the two 
prefactors to the power-law decay is not observed in the 
expansion from a TL state. There, a power-law with a 
single pair of exponent and prefactor governs the decay 
of one-particle and density-density correlations [16]. 



E. Nonintegrable systems 

Beyond the case of the Hubbard model Eq. (1), the 
question arises whether non-stationary states of other 
model Hamiltonians may as well be described by ground- 
state reference systems. Conceptually, one may wonder 
whether integrability plays a role or not. 

While a full account of these interesting issues is be- 
yond the scope of the present work, we wish to at least 
comment on one additional model, the extended Hub- 
bard model. In addition to the terms given in Eq. (1), 
this model incorporates a nearest-neighbor repulsion Hi : 
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FIG. 6: (Color online) Structure factor for (a) spin and (b) 
density correlations for the extended Hubbard model with 
U = 8t, V = t. The results are for the expansion from an 
initial state with ni n i t = 0.6 and 6 particles. Time-dependent 
data are represented with solid lines and solid symbols, the 
corresponding ground state reference system ones are dis- 
played with open circles. 



constructed according to the prescription of Sec. IV A. 
The results for U = 8t and V = t are collected in Fig. 6, 
and we note that in this case, we consider the expansion 
from a box trap only. 

Our observation is that the reference systems still pro- 
vide a very good approximation to the time-dependent 
correlations. Moreover, the agreement remains to be 
much better for the spin structure factor than for the 
charge structure factor. Future work will have to clarify 
whether the reference system description breaks down as 
V is increased. In conclusion, we find that the reference 
system description does not seem to be crucially depen- 
dent on integrability, and we expect a similar picture to 
emerge in other models. 



H 2 = V 



N-l 

^ nini + i 
1=1 



(7) 



Discussion: Relation to density functional 
theory 



The nearest-neighbor interaction both renders the sys- 
tem nonintegrable and induces additional phases at half- 
filling. For V < U/2, the system is a Mott insulator, 
while a large V drives the system into a charge-density 
wave phase [28]. Here we focus on the numerically less 
demanding case of the expansion from an initial state 
with an incommensurate filling of n^t = 0.6. We post- 
pone the discussion of the MDF to a future publication, 
but rather compute the spin and charge structure factors 
for both the expanding system and the reference systems 



Our main results have shown that starting from a 
Mott-insulator, the temporal evolution of spin and den- 
sity correlation functions are very accurately described 
by the ground-state of reference systems defined at each 
instant of time so that they have the same density distri- 
bution as the time evolving system. Such a description 
deteriorates after the Mott insulating region has totally 
melted. On the other hand, for systems starting with 
densities n m i t < 1, the description with reference sys- 
tems remains valid up to the largest times reached in our 
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simulations. These results explicitcly show that the cor- 
relation functions studied here are functionals of the den- 
sity, a fact that is in accordance with density functional 
theory (DFT) for time-dependent systems [29, 30]. 

DFT considers, both for the ground state and time- 
dependent situations, this kind of Hamiltonian [29, 30]: 

H = i?kin + H Coulomb + #ext , (8) 

i.e., one separates the Hamiltonian into kinetic energy, 
the interaction energy due to Coulomb interactions, and 
an external potential. We should here remark that in our 
case, the time-dependent external potential is discontin- 
uous at time r = 0, as we use H axt = H con t at time r = 
and -ffcxt = for r > 0. Therefore, it does not strictly 
comply with the assumptions needed to prove the Rungc- 
Gross theorem [31], i.e., the external potential needs to 
be analytical around r = 0. 

However, and most importantly, the reference sys- 
tems explicitely provide the required functionals, namely 
the correlation functions in the respective ground-states. 
This is, to our opinion, a rather surprising and nontrivial 
fact that the correlations of a genuinely non-stationary 
state can be quantitatively described by equilibrium sys- 
tems. Moreover, since they are in their ground-state, this 
suggests that a minimum principle is at work here. As 
shown in Sec. IV E, these conclusions are not restricted 
to the pure Hubbard-model, i.e., they are not a conse- 
quence of intcgrability. Therefore, we expect that they 
hold in general. 



V. SUMMARY 

In this work, we have identified several remarkable 
and unexpected properties of fcrmions expanding into an 
empty lattice. These include the emergence of coherence 
as well as well as an accumulation of particles at mo- 
mentum 7r/2 in the expansion of particles coming from 
a MI region. In particular, we have shown that correla- 
tion functions of expanding, interacting fermions can be 
accurately described by equilibrium reference systems in 
their ground state. These results are expected to quali- 
tatively carry over to other models as well, and certainly 
also apply to the case of hard-core bosons. 
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